PatchSizePilot
Emanuele Giacomuzzo
2022-08-03
Overview
Data
Meta-ecosystems
All
Equally sized
MM vs SL
SL vs SL*
Patches
All
Isolated
Alpha Diversity
ds_patches %>%
filter(disturbance == disturbance_input,
metaecosystem == "no") %>%
ungroup() %>%
ggplot(
aes(
x = day,
y = shannon,
group = interaction(day, eco_metaeco_type),
fill = eco_metaeco_type
)
) +
geom_boxplot() +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
labs(
x = "Day",
y = "Shannon Index",
fill = "",
caption = "Vertical grey line: resource flow"
) +
scale_fill_manual(
values = c(colour_large, colour_medium, colour_small),
labels = c("Large isolated",
"Medium isolated",
"Small isolated")
) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
)
Small
Alpha Diversity
ds_patches %>%
filter(disturbance == disturbance_input,
patch_size == "S") %>%
ggplot(
aes(
x = day,
y = shannon,
group = culture_ID,
fill = culture_ID,
color = culture_ID,
linetype = eco_metaeco_type
)
) +
geom_line(stat = "summary", fun = "mean") +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
labs(
x = "Day",
y = "Shannon Index",
linetype = "",
caption = "Vertical grey line: resource flow"
) +
scale_linetype_discrete(labels = c(
"Small isolated",
"Small connected to large",
"Small connected to small"
)) +
scale_x_continuous(breaks = unique(ds_patches$day)) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
)
ds_patches %>%
filter(disturbance == disturbance_input,
patch_size == "S") %>%
ggplot(aes(
x = day,
y = shannon,
group = interaction(day, eco_metaeco_type),
fill = eco_metaeco_type
)) +
geom_boxplot() +
labs(
x = "Day",
y = "Alpha Diversity (Shannon Index)",
fill = "",
caption = "Vertical grey line: resource flow"
) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
) +
scale_fill_manual(
values = c(colour_isolated, colour_different_size, colour_same_size),
labels = c(
"Small isolated",
"Small connected to large",
"Small connected to small"
)
) +
scale_x_continuous(breaks = unique(ds_patches$day)) +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
)
ds_patches_effect_size %>%
filter(
!time_point == 0,
disturbance == disturbance_input,
eco_metaeco_type %in% c("S (S_S)", "S (S_L)")
) %>%
ggplot(aes(x = day,
y = shannon_d,
color = eco_metaeco_type)) +
geom_point(position = position_dodge(0.5)) +
geom_line(position = position_dodge(0.5)) +
geom_errorbar(
aes(ymin = shannon_d_lower,
ymax = shannon_d_upper),
width = .2,
position = position_dodge(0.5)
) +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_hline(
yintercept = 0,
linetype = "dotted",
color = "black",
size = 0.7
) +
labs(
title = paste("Disturbance =", disturbance_input),
x = "Day",
y = "Shannon Index Effect Size (Hedge's d)",
color = ""
) +
scale_color_discrete(labels = c("Small connected to large",
"Small connnected to small")) +
scale_x_continuous(breaks = unique(ds_patches_effect_size$day)) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
) 
Let’s test our hypothesis that small patches connected to large had more species richness. To do so, we will use a linear model. As we are using an effect size (control = isolated small), we have pulled all the replicates together and we cannot include random effects. We will exclude time point 0 because there all the cultures were the same (we did not measure each culture at day 0, but we measured the large bottle from which we created all cultures). We excluded the data points of time point = 1 because disturbance and ecosystem type couldn’t have had an effect, but we included them as baselines. We are now including both baselines and their interaction with time, as baselines should have less of an effect as time goes by.
\[ Species \: Richness \: (Hedge's \: d) = B + t + P + D + t*P + t*D + P*D + t*B \]
B = baseline_shannon (bioarea at t1, before the first resource flow)
t = time
P = patch type (connected to small vs connected to large)
D = disturbance
filtered_data = ds_patches_effect_size %>%
filter(time_point >= first_time_point_model,
time_point <= last_time_point_model,
disturbance == disturbance_input,
eco_metaeco_type %in% c("S (S_S)", "S (S_L)"))full_model = lm(shannon_d ~
day +
eco_metaeco_type +
day : eco_metaeco_type +
baseline_shannon +
day : baseline_shannon,
data = filtered_data)
null_model = lm(shannon_d ~
day,
data = filtered_data)
anova(full_model, null_model)## Analysis of Variance Table
##
## Model 1: shannon_d ~ day + eco_metaeco_type + day:eco_metaeco_type + baseline_shannon +
## day:baseline_shannon
## Model 2: shannon_d ~ day
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 6 4.8394
## 2 10 16.5368 -4 -11.697 3.6257 0.07824 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(full_model, null_model)## df AIC
## full_model 7 37.15707
## null_model 3 43.90272Let’s then do some model diagnostics.
par(mfrow = c(2,3))
plot(full_model, which = 1:5)
R2_full = glance(full_model)$r.squared
R2_null = glance(null_model)$r.squared
R2_P = R2_full - R2_null
R2_full = round(R2_full, digits = 2)
R2_P = round(R2_P, digits = 2)The adjusted R squared of the model is 0.78 and the adjusted R squared of patch type is 0.53 (which includes also its interaction with disturbance). I don’t think it’s necessary to look at effect sizes at different time points.
Large
Alpha Diversity
ds_patches %>%
filter(disturbance == disturbance_input,
patch_size == "L") %>%
ggplot(
aes(
x = day,
y = shannon,
group = culture_ID,
fill = culture_ID,
color = culture_ID,
linetype = eco_metaeco_type
)
) +
geom_line(stat = "summary", fun = "mean") +
labs(
x = "Day",
y = "Shannon Index",
linetype = "",
caption = "Vertical grey line: resource flow"
) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
) +
scale_linetype_discrete(labels = c(
"Large isolated",
"Large connected to large",
"Large connected to small"
)) +
scale_x_continuous(breaks = unique(ds_patches$day)) +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
)
ds_patches %>%
filter(disturbance == disturbance_input,
patch_size == "L") %>%
ggplot(aes(
x = day,
y = shannon,
group = interaction(day, eco_metaeco_type),
fill = eco_metaeco_type
)) +
geom_boxplot() +
labs(
title = paste("Disturbance =", disturbance_input),
x = "Day",
y = "Shannon Index",
fill = "",
caption = "Vertical grey line: resource flow"
) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
) +
scale_fill_manual(
values = c(colour_isolated, colour_same_size, colour_different_size),
labels = c(
"Small isolated",
"Large connected to large",
"Large connected to small"
)
) +
scale_x_continuous(breaks = unique(ds_patches$day)) +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
)
ds_patches_effect_size %>%
filter(
!time_point == 0,
disturbance == disturbance_input,
eco_metaeco_type == "L (L_L)" |
eco_metaeco_type == "L (S_L)"
) %>%
ggplot(aes(x = day,
y = shannon_d,
color = eco_metaeco_type)) +
geom_point(position = position_dodge(0.5)) +
geom_line(position = position_dodge(0.5)) +
geom_errorbar(
aes(ymin = shannon_d_lower,
ymax = shannon_d_upper),
width = .2,
position = position_dodge(0.5)
) +
geom_vline(
xintercept = resource_flow_days[1] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[2] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[3] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[4] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[5] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_vline(
xintercept = resource_flow_days[6] + 0.7,
linetype = "dotdash",
color = "grey",
size = 0.7
) +
geom_hline(
yintercept = 0,
linetype = "dotted",
color = "black",
size = 0.7
) +
scale_x_continuous(breaks = unique(ds_patches_effect_size$day)) +
labs(
title = paste("Disturbance =", disturbance_input),
x = "Day",
y = "Shannon Index Effect Size (Hedge's d)",
color = ""
) +
scale_color_discrete(labels = c("Large connected to large",
"Large connected to small")) +
theme_bw() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "bottom"
)
Let’s test our hypothesis that small patches connected to large had more species richness. To do so, we will use a linear model. As we are using an effect size (control = isolated small), we have pulled all the replicates together and we cannot include random effects. We will exclude time point 0 because there all the cultures were the same (we did not measure each culture at day 0, but we measured the large bottle from which we created all cultures). We excluded the data points of time point = 1 because disturbance and ecosystem type couldn’t have had an effect, but we included them as baselines. We are now including both baselines and their interaction with time, as baselines should have less of an effect as time goes by.
\[ Species \: Richness \: (Hedge's \: d) = B + t + P + D + t*P + t*D + P*D + t*B \]
B = baseline_shannon (bioarea at t1, before the first resource flow)
t = time
P = patch type (connected to small vs connected to large)
D = disturbance
filtered_data = ds_patches_effect_size %>%
filter(time_point >= first_time_point_model,
time_point <= last_time_point_model,
disturbance == disturbance_input,
eco_metaeco_type %in% c("S (S_S)", "S (S_L)"))full_model = lm(shannon_d ~
day +
eco_metaeco_type +
day : eco_metaeco_type,
data = filtered_data)
null_model = lm(shannon_d ~
day,
data = filtered_data)
anova(full_model, null_model)## Analysis of Variance Table
##
## Model 1: shannon_d ~ day + eco_metaeco_type + day:eco_metaeco_type
## Model 2: shannon_d ~ day
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 8 11.232
## 2 10 16.537 -2 -5.305 1.8893 0.2128AIC(full_model, null_model)## df AIC
## full_model 5 43.26070
## null_model 3 43.90272Let’s then do some model diagnostics.
par(mfrow = c(2,3))
plot(full_model, which = 1:5)
R2_full = glance(full_model)$r.squared
R2_null = glance(null_model)$r.squared
R2_P = R2_full - R2_null
R2_full = round(R2_full, digits = 2)
R2_P = round(R2_P, digits = 2)The adjusted R squared of the model is 0.49 and the adjusted R squared of patch type is 0.24 (which includes also its interaction with disturbance). I don’t think it’s necessary to look at effect sizes at different time points.
Evaporation
Effects of evaporation
We want to know if there was a systematic bias in the evaporation of different treatments (disturbance, patch size) and whether evaporation changed across time. My expectation would be that we would see a difference among the exchanges 2,3 and the exchanges 4,5,6. This is because in exchange 2,3 cultures were microwaved in 15 tubes for 3 minutes and in exchange 4,5,6 cultures were microwaved in 4 tubes for only 1 minute.
Tidy
#Columns: exchange & evaporation
ds_for_evaporation = gather(ds_for_evaporation,
key = exchange,
value = evaporation,
water_add_after_t2:water_add_after_t6)
ds_for_evaporation[ds_for_evaporation == "water_add_after_t2"] = "2"
ds_for_evaporation[ds_for_evaporation == "water_add_after_t3"] = "3"
ds_for_evaporation[ds_for_evaporation == "water_add_after_t4"] = "4"
ds_for_evaporation[ds_for_evaporation == "water_add_after_t5"] = "5"
ds_for_evaporation[ds_for_evaporation == "water_add_after_t6"] = "6"
ds_for_evaporation$evaporation[ds_for_evaporation$exchange == 2] = ds_for_evaporation$evaporation[ds_for_evaporation$exchange == 2] / 2 #This is because exchange contained the topping up of two exchanges
ds_for_evaporation$evaporation[ds_for_evaporation$exchange == 2] = ds_for_evaporation$evaporation[ds_for_evaporation$exchange == 2] + 2 #We need to add 2 ml to the evaporation that happened at the exchange events 1 and 2. This is because we already added 1 ml of water at exchange 1 and 1 ml of water at exchange 2.
#Column: nr_of_tubes_in_rack
ds_for_evaporation$nr_of_tubes_in_rack[ds_for_evaporation$exchange == 1] = 15
ds_for_evaporation$nr_of_tubes_in_rack[ds_for_evaporation$exchange == 2] = 15
ds_for_evaporation$nr_of_tubes_in_rack[ds_for_evaporation$exchange == 3] = 15
ds_for_evaporation$nr_of_tubes_in_rack[ds_for_evaporation$exchange == 4] = 4
ds_for_evaporation$nr_of_tubes_in_rack[ds_for_evaporation$exchange == 5] = 4
ds_for_evaporation$nr_of_tubes_in_rack[ds_for_evaporation$exchange == 6] = 4Plot
ds_for_evaporation %>%
filter(disturbance == disturbance) %>%
ggplot(aes(x = as.character(nr_of_tubes_in_rack),
y = evaporation)) +
geom_boxplot() +
labs(x = "Number of tubes in rack",
y = "Evaporation (ml)")
ds_for_evaporation %>%
filter(disturbance == disturbance) %>%
ggplot(aes(x = as.character(patch_size),
y = evaporation)) +
geom_boxplot() +
labs(x = "Patch size",
y = "Evaporation (ml)")
ds_for_evaporation %>%
filter(disturbance == disturbance) %>%
ggplot(aes(x = as.character(day),
y = evaporation)) +
geom_boxplot() +
labs(x = "Day",
y = "Evaporation (ml)")
ds_for_evaporation %>%
filter(disturbance == disturbance) %>%
ggplot(aes(x = disturbance,
y = evaporation)) +
geom_boxplot() +
labs(x = "Disturbance",
y = "Evaporation (ml)")
It seems like there is no real difference across time, disturbance, or patch type. However, we could also run a mixed effect model to show that they do not.
Mixed effect model
It gives me the following error:
- Error in fn(nM$xeval()) : Downdated VtV is not positive definite
mixed.model = lmer(evaporation ~
patch_size * disturbance * exchange +
(exchange | culture_ID),
data = ds_for_evaporation,
REML = FALSE,
control = lmerControl (optimizer = "Nelder_Mead"))
null.model = lm(evaporation ~
1,
data = ds_for_evaporation)
anova(mixed.model, null.model)Evaporation tests
Evaporation when microwaving 15 falcon tubes at the time
evaporation.test = read.csv(here("data", "evaporation_test","evaporation_test_right.csv"), header = TRUE)
evaporation.test %>%
ggplot(aes (x = as.character(water_pipetted),
y = weight_water_evaporated,
group = interaction(water_pipetted, as.character(rack)),
fill = as.character(rack))) +
geom_boxplot() +
labs(x = "Water volume (ml)" ,
y = "Evaporation (g)",
fill = "Rack replicate")
Evaporation when microwaving 5 tubes with 10 filled or empty tubes
evaporation.test = read.csv(here("data", "evaporation_test", "evaporation_test_fill_nofill.csv"), header = TRUE)
evaporation.test %>%
ggplot(aes (x = all_tubes_water,
y = weight_water_evaporated)) +
geom_boxplot() +
labs(x = "Water in the other 10 tubes" ,
y = "Evaporation (g)",
caption = "When all tubes were filled, they were filled with 6.75 ml of deionised water.")
Video analysis script
rm(list = ls())
setwd("/media/mendel-himself/ID_061_Ema2/PatchSizePilot/training")
# load package
# library(devtools)
# install_github("femoerman/bemovi", ref="master")
library(bemovi)
library(parallel)
library(doParallel)
library(foreach)
#Define memory to be allocated
memory.alloc <- 240000 #-needs_to_be_specified
memory.per.identifier <- 40000 #-needs_to_be_specified
memory.per.linker <- 5000 #-needs_to_be_specified
memory.per.overlay <- 60000 #-needs_to_be_specified
# UNIX
# set paths to tools folder and particle linker
tools.path <-
"/home/mendel-himself/bemovi_tools/" #-needs_to_be_specified
to.particlelinker <- tools.path
# directories and file names
to.data <- paste(getwd(), "/", sep = "")
video.description.folder <- "0_video_description/"
video.description.file <- "video_description.txt"
raw.video.folder <- "1_raw/"
raw.avi.folder <- "1a_raw_avi/"
metadata.folder <- "1b_raw_meta/"
particle.data.folder <- "2_particle_data/"
trajectory.data.folder <- "3_trajectory_data/"
temp.overlay.folder <- "4a_temp_overlays/"
overlay.folder <- "4_overlays/"
merged.data.folder <- "5_merged_data/"
ijmacs.folder <- "ijmacs/"
######################################################################
# VIDEO PARAMETERS
# video frame rate (in frames per second)
fps <- 25 #-needs_to_be_specified
# length of video (in frames)
total_frames <- 125 #-needs_to_be_specified
#Dimensions of the videos in pixels
width = 2048 #-needs_to_be_specified
height = 2048 #-needs_to_be_specified
# measured volume (in microliter) #-needs_to_be_specified
measured_volume <-
34.4 # for Leica M205 C with 1.6 fold magnification, sample height 0.5 mm and Hamamatsu Orca Flash 4
#measured_volume <- 14.9 # for Nikon SMZ1500 with 2 fold magnification, sample height 0.5 mm and Canon 5D Mark III
# size of a pixel (in micrometer) #-needs_to_be_specified
pixel_to_scale <-
4.05 # for Leica M205 C with 1.6 fold magnification, sample height 0.5 mm and Hamamatsu Orca Flash 4
#pixel_to_scale <- 3.79 # for Nikon SMZ1500 with 2 fold magnification, sample height 0.5 mm and Canon 5D Mark III
# specify video file format (one of "avi","cxd","mov","tiff")
# bemovi only works with avi and cxd. other formats are reformated to avi below
video.format <- "cxd" #-needs_to_be_specified
# setup
difference.lag <- 10
thresholds <- c(13, 255) # don't change the second value
# thresholds <- c(50,255)
# MORE PARAMETERS (USUALLY NOT CHANGED)
######################################################################
# FILTERING PARAMETERS
# optimized for Perfex Pro 10 stereomicrocope with Perfex SC38800 (IDS UI-3880LE-M-GL) camera
# tested stereomicroscopes: Perfex Pro 10, Nikon SMZ1500, Leica M205 C
# tested cameras: Perfex SC38800, Canon 5D Mark III, Hamamatsu Orca Flash 4
# tested species: Tet, Col, Pau, Pca, Eug, Chi, Ble, Ceph, Lox, Spi
# min and max size: area in pixels
particle_min_size <- 10
particle_max_size <- 1000
# number of adjacent frames to be considered for linking particles
trajectory_link_range <- 3
# maximum distance a particle can move between two frames
trajectory_displacement <- 16
# these values are in the units defined by the parameters above: fps (seconds), measured_volume (microliters) and pixel_to_scale (micometers)
filter_min_net_disp <- 25
filter_min_duration <- 1
filter_detection_freq <- 0.1
filter_median_step_length <- 3
######################################################################
# VIDEO ANALYSIS
#Check if all tools are installed, and if not install them
check_tools_folder(tools.path)
#Ensure computer has permission to run bftools
system(paste0("chmod a+x ", tools.path, "bftools/bf.sh"))
system(paste0("chmod a+x ", tools.path, "bftools/bfconvert"))
system(paste0("chmod a+x ", tools.path, "bftools/showinf"))
# Convert files to compressed avi (takes approx. 2.25 minutes per video)
convert_to_avi(
to.data,
raw.video.folder,
raw.avi.folder,
metadata.folder,
tools.path,
fps,
video.format
)
# TESTING
# check file format and naming
# check_video_file_names(to.data,raw.avi.folder,video.description.folder,video.description.file)
# check whether the thresholds make sense (set "dark backgroud" and "red")
#check_threshold_values(to.data, raw.avi.folder, ijmacs.folder, 2, difference.lag, thresholds, tools.path, memory.alloc)
# identify particles
locate_and_measure_particles(
to.data,
raw.avi.folder,
particle.data.folder,
difference.lag,
min_size = particle_min_size,
max_size = particle_max_size,
thresholds = thresholds,
tools.path,
memory = memory.alloc,
memory.per.identifier = memory.per.identifier,
max.cores = detectCores() - 1
)
# link the particles
link_particles(
to.data,
particle.data.folder,
trajectory.data.folder,
linkrange = trajectory_link_range,
disp = trajectory_displacement,
start_vid = 1,
memory = memory.alloc,
memory_per_linkerProcess = memory.per.linker,
raw.avi.folder,
max.cores = detectCores() - 1,
max_time = 1
)
# merge info from description file and data
merge_data(
to.data,
particle.data.folder,
trajectory.data.folder,
video.description.folder,
video.description.file,
merged.data.folder
)
# load the merged data
load(paste0(to.data, merged.data.folder, "Master.RData"))
# filter data: minimum net displacement, their duration, the detection frequency and the median step length
trajectory.data.filtered <-
filter_data(
trajectory.data,
filter_min_net_disp,
filter_min_duration,
filter_detection_freq,
filter_median_step_length
)
# summarize trajectory data to individual-based data
morph_mvt <-
summarize_trajectories(
trajectory.data.filtered,
calculate.median = F,
write = T,
to.data,
merged.data.folder
)
# get sample level info
summarize_populations(
trajectory.data.filtered,
morph_mvt,
write = T,
to.data,
merged.data.folder,
video.description.folder,
video.description.file,
total_frames
)
# create overlays for validation
create.subtitle.overlays(
to.data,
traj.data = trajectory.data.filtered,
raw.video.folder,
raw.avi.folder,
temp.overlay.folder,
overlay.folder,
fps,
vid.length = total_frames / fps,
width,
height,
tools.path = tools.path,
overlay.type = "number",
video.format
)
# Create overlays (old method)
create_overlays(
traj.data = trajectory.data.filtered,
to.data = to.data,
merged.data.folder = merged.data.folder,
raw.video.folder = raw.avi.folder,
temp.overlay.folder = "4a_temp_overlays_old/",
overlay.folder = "4_overlays_old/",
width = width,
height = height,
difference.lag = difference.lag,
type = "traj",
predict_spec = F,
contrast.enhancement = 1,
IJ.path = "/home/mendel-himself/bemovi_tools",
memory = memory.alloc,
max.cores = detectCores() - 1,
memory.per.overlay = memory.per.overlay
)
########################################################################
# some cleaning up
#system("rm -r 2_particle_data")
#system("rm -r 3_trajectory_data")
#system("rm -r 4a_temp_overlays")
system("rm -r ijmacs")
########################################################################Here I want to check whether I can distinguish between Tet and Chi (food of Ble). To do so, let’s look at the traits of Ble and Chi.
load(here("data", "morphology", "training.RData"))
training = morph_mvt
rm(morph_mvt)training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_grey)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_area)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_perimeter)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_major)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_minor)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_ar)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = mean_turning)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = duration)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = N_frames)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = max_net)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = net_disp)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = net_speed)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = gross_disp)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = gross_speed)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = max_step)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = min_step)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = max_gross_speed)) +
geom_boxplot()
training %>%
filter(species == "Ble" | species == "Tet") %>%
ggplot(aes(x = species,
y = min_gross_speed)) +
geom_boxplot()
It seems like there’s no great way of telling them apart.
Running time
## Time difference of 1.426613 mins